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ABSTRACT 

This  project  focused  on  the  development  of  multiscale  models  for  molecular  crystals,  with  special  emphasis  on  the  secondary  explosive 
cyclotrimethylene  trinitramine  (RDX).  Physical  understanding  related  to  the  molecular  scale  mechanics  of  this  material  was  also  developed. 
This  data  is  used  as  starting  point  for  the  definition  of  physically-based  coarse  grained  models.  In  turn,  these  models  are  to  assist  the  design 
of  new  energetic  materials  with  enhanced  energy  release  rates  and  reduced  sensitivity  to  unintentional  detonation.  The  following  results  have 
been  obtained:  a)  we  investigated  conformational  stability  of  the  RDX  molecule  in  the  crystal  as  a  function  of  temperature  and  determined 
the  density  of  conformers  in  steady  state  at  temperatures  between  room  and  melting,  b)  we  identified  the  stable  dislocations  in  RDX  and,  for 
the  first  time,  established  a  ranking  of  slip  systems  in  terms  of  the  Peierls  stress  (critical  stress  for  dislocation  motion),  c)  we  discovered  a 
family  of  point  defects  which  are  rotated  and  distorted  molecules  embedded  in  an  otherwise  perfect  crystal;  the  stability  of  such  point  defects 
was  demonstrated  with  DFT  simulations  and  a  method  was  proposed  to  detect  them  experimentally  using  Raman  spectroscopy,  d)  we  are 
currently  investigating  the  effect  of  temperature  on  the  ranking  of  slip  systems,  which  is  a  result  needed  in  crystal  plasticity  models  of  hot 
spot  formation,  e)  we  developed  a  family  of  coarse  grained  models  for  the  RDX  crystal  in  which  the  molecules  are  gradually  coarse  grained 
to  higher  and  higher  degrees,  ranging  from  full  atomistic  detail  to  rigid  blobs,  and  we  determined  the  level  of  error  introduced  by  coarse 
graining  in  various  macroscopic  measures  of  crystal  thermo-mechanics. 
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Scientific  Progress 


The  main  results  of  this  work  are: 


a)  identified  conformational  stability  of  RDX  molecules  in  the  crystal  and  its  dependence  on  temperature, 

b)  identified  the  stable  dislocations  in  RDX, 

c)  established  a  ranking  of  slip  systems  in  terms  of  the  Peierls  stress  at  zero  Kelvin.  We  are  currently  investigating  the  effect  of 
temperature  on  this  ranking, 

d)  discovered  a  family  of  rotational  defects  (point  defects)  in  RDX, 

e)  developed  a  family  of  coarse  grained  models  of  the  RDX  crystal. 

We  outline  briefly  the  central  results  and  their  significance.  The  discussion  here  refers  to  the  attached  file  which  contains  actual 
data  and  schematic  representations. 

a)  Conformational  stability  of  RDX  molecules 

The  RDX  molecule  exists  in  the  form  of  4  conformers  denoted  by  Caae,  Caee,  Caaa  and  Ceee.  The  molecule  is  composed 
from  a  triazine  ring  and  three  nitro  side  groups  (Slide  1).  The  side  groups  take  positions  that  make  angles  ranging  from 
approximately  35o  to  135o  with  respect  to  the  normal  of  the  ring.  An  angle  of  90o  corresponds  to  an  equatorial  position 
(denoted  by  ‘e’),  while  angles  close  to  the  extremes  of  the  indicated  range  correspond  to  axial  positions  (denoted  by  ‘a’).  A 
conformer  of  Caae  type  has  one  nitro  group  in  the  equatorial  position,  and  two  nitro  groups  in  the  axial  position  (Slide  2).  The 
transition  of  a  nitro  group  between  a  and  e  states  is  limited  by  potential  barriers  (Slide  2). 

Relevance  of  conformer  transitions  in  the  crystal.  We  determined  that  the  molecules  do  not  change  state  in  the  crystal  at  room 
temperature  and  in  absence  of  loads.  Transitions  take  place  when  either  tension  (uniaxial  or  biaxial)  or  compression  are 
applied.  Under  such  conditions,  molecules  change  state  by  switching  from  Caae  (the  ground  state)  to  Caee  and  back.  In  doing 
so,  they  absorb  phonons  with  energy  at  least  equal  to  the  height  of  the  barrier  between  the  two  states.  This  energy  is  released 
in  the  form  of  molecular  vibrations  when  the  structure  “falls  from  the  height  of  the  barrier”  into  the  new  potential  well 
(corresponding  to  the  conformer).  The  same  process  happens  when  the  molecule  transforms  back  from  Caee  to  the  ground 
state.  We  have  seen  (Slide  5)  that  this  leads  to  localized  heating  of  the  respective  molecule  by  ~10K.  Hence,  the  molecule 
temperature  increases  immediately  after  the  transition  and  then  decreases  gradually  (ps  time  scale)  by  thermalization. 

Effect  of  strain/stress.  It  has  been  observed  that  conformers  are  not  stable  in  the  unloaded  alpha-RDX  crystal.  However,  we 
determined  that  a  fraction  of  the  Caae  molecules  transform  to  Caee  if  the  crystal  is  stretched  in  the  [001]  crystallographic 
direction.  Interestingly,  stretching  in  other  crystal  directions,  as  well  as  applying  a  hydrostatic  strain  does  not  promote  the 
transformation  (Slide  3). 

Correlated  dynamics  of  conformal  transitions.  Since  a  conformer  induces  a  long-range  strain  field,  neighboring  conformers 
interact.  We  evaluated  the  range  of  interaction  (which  is  quite  large:  conformers  interact  even  when  separated  by  17A)  and 
determined  that  a)  conformers  cluster  in  space  (i.e.  transition  from  the  ground  state  to  Caee  is  favorable  in  the  vicinity  of  a 
molecule  which  has  already  transformed),  b)  conformational  changes  between  Caae  and  Caee  are  correlated  in  time  (Slide  4). 
Both  spatial  and  temporal  correlations  are  attributed  to  the  field-mediated  interaction  of  these  “point  defects.”  Observation  a)  is 
relevant  for  the  concentration  of  the  kinetic  (thermal)  energy  on  the  molecular  scale.  As  discussed  above,  when  an  isolated 
molecule  transforms  between  conformer  states,  its  effective  temperature  (kinetic  energy)  increases  by  approximately  10K.  If 
multiple  molecules  in  a  region  of  the  crystal  perform  similar  transitions,  the  temperature  of  a  larger  region  of  the  crystal 
increases.  This  fluctuation  lasts  longer  than  the  temperature  fluctuation  resulting  when  a  single  molecule  transforms  since 
thermalization  takes  longer  when  a  larger  volume  of  material  is  involved. 

b)  Stability  of  dislocations  in  RDX 

Initiation  and  detonation  happen  as  a  consequence  of  plastic  deformation  and/or  overheating.  In  this  work  we  investigated  the 
mechanisms  by  which  plastic  deformation  takes  place  in  RDX. 

RDX,  with  its  complex  crystal  and  molecular  structure,  is  generally  considered  elastically  soft,  but  plastically  hard.  This  is  a 
consequence  of  the  relatively  soft  inter-molecular  bonding  and  is  a  common  feature  in  molecular  crystals.  It  has  been  shown 
experimentally  that  the  easiness  of  slip  in  certain  directions  correlates  with  initiation.  In  this  work  we  aimed  to  identify  the  active 
slip  systems  and  to  rank  them  in  terms  of  dislocation  mobility  such  to  help  make  this  statement  quantitative. 

It  results  that  in  many  slip  systems  in  which  dislocations  could  -in  principle-  form,  plastic  deformation  is  actually  sluggish  since 
either  only  dislocations  with  large  Burgers  vectors  exist,  or  the  stable  dislocations  have  high  Peierls  stresses.  To  clarify  on 
which  crystallographic  planes  and  directions  dislocations  can  form,  we  focused  on  investigating  the  core-structures  of  edge  and 
screw  dislocations  in  various  experimentally  reported  slip  systems  for  RDX.  Experimental  studies  suggest  the  occurrence  of  slip 
on  the  {021},  {011}  and  (010)  planes.  In  addition  to  these,  the  slip  plane  (001)  was  also  considered.  The  crystal  structure  of 
RDX  permits  slip  in  two  parallel  planes  for  each  type.  In  our  investigations,  the  plane  with  lower  attachment  energy  was 
selected.  Calculations  of  the  dislocation  energies  were  done  using  an  atomistic  model  using  the  Smith-Bharadwaj  potential. 

Slide  6  shows  a  (100)  projection  of  the  RDX  crystal  in  its  ground  state  and  the  slip  systems  which  produce  slip  in  the  [100] 
direction.  Note  that  several  glide  planes  may  produce  this  type  of  slip  and  the  associated  dislocations  can  be  edge  or  screw. 
Slides  7,  8  and  9  show  core  structures  for  various  glide  planes  and  dislocation  types.  The  cores  are  generally  narrow  and 
splitting  of  dislocations  in  partials  is  rarely  observed.  The  narrowest  core  is  seen  for  the  (01 1 )  [1 00]  system  with  a  width  of 


approximately  3*1100.,  where  I  stands  for  the  unit  cell  in  the  respective  direction.  The  (010)  and  (021)  edge  dislocations  in  the 
[100]  direction  have  fairly  wide  cores  with  a  width  of  6*  1100.  An  interesting  feature  can  be  observed  for  the  (001 )  [010]  edge 
dislocation  where  the  stable  core  has  a  width  of  about  5*1010  and  develops  out  of  plane  screw  component. 

c)  Ranking  of  slip  systems 

The  central  result  related  to  dislocation  mechanics,  which  is  the  physical  basis  of  plasticity,  is  the  ranking  of  slip  systems.  The 
potentially  active  slip  systems  have  been  identified  using  atomistic  models  (section  (b)),  but  some  of  these  systems  host 
essentially  sessile  dislocations  and  are  effectively  inactive.  To  clarify  which  slip  systems  are  actually  active  one  has  to  evaluate 
the  Peierls  stress,  i.e.  the  stress  required  to  move  dislocations  in  the  perfect  crystal.  These  stress  values  represent  the  lower 
limit  of  the  stress  required  for  plastic  deformation  in  the  various  crystallographic  planes  and  directions.  They  can  be  used  to 
determine  which  systems  are  active  and  whether  a  sufficient  number  of  systems  exist  to  render  plastic  deformation  feasible 
macroscopically.  The  importance  of  this  result  from  a  modeling/computation  standpoint  is  that  in  coarse  grained  dislocation 
dynamics  (DDD)  simulations,  one  needs  to  input  the  Peierls  stress  which  enters  the  equation  of  motion  of  dislocations.  Hence, 
the  study  performed  here  represents  a  sequential  scale  linking  step  in  this  multiscale  analysis  of  plasticity  in  RDX. 

Slide  10  shows  the  Peierls  stress  values  calculated  using  a  large  scale  atomistic  model.  The  4  types  of  symbols  correspond  to  4 
slip  planes  which  are  considered  to  be  active  in  RDX.  In  each  plane  one  may  have  edge  and  screw  dislocations  and  these  may 
have  full  or  half  Burger  vector.  The  notation  “e1/2[100]”,  for  example,  stands  for  the  edge  dislocation  with  Burgers  vector  14  of 
the  full  Burgers  vector  in  the  [100]  direction  of  the  crystal.  The  horizontal  axis  shows  the  theoretical  strength  of  the  respective 
systems  computed  based  on  the  shear  modulus  in  the  respective  crystal  direction.  It  is  seen  that  the  system  with  lowest  Peierls 
barrier  is  the  14[100]  system,  in  which  both  the  screw  and  edge  dislocations  have  critical  stresses  about  an  order  of  magnitude 
smaller  than  the  theoretical  strength. 

Slide  1 1  shows  some  of  the  implications  of  the  ranking  of  slip  systems  in  Slide  10.  The  figure  on  the  left  shows  the  resolved 
stress  in  all  systems  activated  by  shocks  applied  in  specific  crystal  directions  and  of  specific  magnitudes.  Shock  of  three 
magnitudes  (from  1 .24  to  7  GPa)  have  been  applied  in  experiments  perpendicular  to  the  (111)  plane.  The  symbols  indicate  the 
resolved  stress  versus  the  Peierls  stress  for  all  systems  activated  by  this  shock.  The  line  has  slope  1 ,  such  that  if  a  point  is 
above  the  line,  the  respective  system  is  active.  It  is  seen  that  for  a  shock  of  magnitude  1 .24  GPa,  all  points  are  below  the  line 
and  hence  there  is  no  plasticity.  For  a  shock  of  7  GPa  all  points  are  above  the  line  and  one  has  plastic  deformation  in  all 
systems.  In  experiments  it  is  seen  that  detonation  happens  when  the  shock  magnitude  is  larger  than  2  GPa.  Here  we  see  that 
when  this  parameter  increases  from  1.24  to  2.14  GPa,  the  points  move  above  the  line  and  more  systems  are  activated.  Hence, 
the  data  are  in  good  agreement  with  the  experimental,  macroscopic  observations  and  suggest  that  plasticity  is  intimately  related 
to  detonation.  The  figure  on  the  right  shows  a  similar  comparison,  with  a  set  of  indentation  experiments  published  in  the 
literature.  The  three  types  of  symbols  correspond  to  indentations  in  different  directions.  The  data  points  correspond  to  slip 
systems  that  may  be  activated  by  the  respective  indentation.  If  the  point  is  above  the  line,  the  respective  system  is  indeed 
activated,  i.e.  the  applied  stress  is  larger  than  the  Peierls  stress.  These  results  are  compared  with  direct  observations  of 
activation  of  plastic  deformation  in  specific  clip  systems  in  these  experiments. 

We  are  in  the  process  of  evaluating  the  temperature  dependence  of  this  ranking.  The  work  will  be  completed  after  the  end  date 
of  the  current  project  and  is  part  of  the  PhD  thesis  of  one  of  the  students  supported  from  the  project  (Anirban  Pal).  It  may  turn 
out  that  as  the  temperature  increases,  some  slip  systems  which  are  not  active  at  temperatures  close  to  zero  Kelvin,  become 
more  active  than  others  that  proved  to  have  lower  Peierls  stress  at  OK.  Hence,  due  to  the  fact  that  the  energetic  barriers 
involved  have  different  height,  the  ranking  at  room  temperature  may  be  quite  different  from  that  reported  here  (Slide  10).  This 
information  is  central  for  the  calibration  of  crystal  plasticity  models. 

Slip  asymmetry 

Slip  asymmetry  is  a  common  occurrence  in  some  monatomic  crystals  where  it  is  due  to  complex  core  structures  or  specific 
packing  of  slip  planes.  We  identified  a  new  mechanism,  based  on  molecular  steric  hindrance,  which  leads  to  asymmetric 
dislocation  motion  in  RDX.  Dislocations  move  at  different  critical  stresses  when  shear  is  applied  in  the  positive  and  negative 
directions  of  the  Burgers  vector  in  the  slip  system  that  contributes  most  to  plastic  deformation.  This  is  a  new  physical 
observation  which  has  significant  implications  for  coarse  grained  DDD  modeling  of  plastic  deformation. 

Slide  12  shows  a  close-up  view  of  the  crystal  in  the  [001]  projection.  The  glide  planes  of  interest  are  horizontal  in  this  view. 
Molecules  are  shown  by  lines  connecting  the  atomic  positions,  while  spheres  represent  the  positions  of  the  center  of  mass  of 
the  molecules.  It  is  clear  from  the  COM  representation  that  there  are  two  distinct  slip  planes  of  the  (010)  type,  which  differ  in 
inter-molecular  bonding  across  the  glide  plane.  These  planes  are  denoted  as  PI  and  P2.  The  interface  between  the  red  and 
blue  layers  is  the  PI  slip  plane  and  the  interface  between  the  orange  and  blue  layers  is  the  P2  plane.  Across  the  PI  plane  the 
nitro  group  of  the  molecules  below  is  oriented  in  the  first  trigonometric  quadrant.  However,  it  is  observed  that  the  PI  plane 
below  the  P2  plane  indicated  in  the  figure  is  not  identical  to  the  PI  plane  above  P2.  In  this  case,  the  nitro  groups  are  pointing  in 
the  second  quadrant.  Therefore,  one  has  two  types  of  PI  planes  which  are  denoted  as  PI  and  PT.  A  similar  situation  exists  for 
the  other  plane  and  two  variants,  P2  and  P2’,  are  defined. 

The  simple  inspection  of  the  crystal  geometry  is  sufficient  to  suggest  the  existence  of  slip  asymmetry.  To  quantify  the  effect,  we 
compute  the  Peierls  stress  for  dislocations  moving  in  the  two  directions  for  all  these  (010)  planes.  The  [100]  direction  pointing  to 
the  right  in  Fig.  (b)  on  slide  12,  is  denoted  as  S+  and  the  opposite  direction  is  denoted  as  S-.  Motion  in  S+  requires  moving  the 


upper  molecular  planes  to  the  right  relative  to  the  lower  planes,  in  this  figure.  The  resolved  shear  stress  required  for  motion  in 
the  S+  and  S-  directions  in  the  PI  plane  is  0.157  GPa  and  0.332  GPa,  respectively,  with  an  uncertainty  of  0.017  GPa.  This 
indicates  that  dislocation  motion  in  these  planes  is  asymmetric  if  the  applied  stress  level  is  below  the  largest  values  of  the 
Peierls  stress.  The  similar  analysis  performed  for  the  two  directions  of  dislocation  motion  in  the  P2  plane  leads  to  a  smaller 
Peierls  stress  and  a  smaller  asymmetry.  Specifically,  we  obtain  0.201  GPa  and  0.263GPa  with  an  uncertainty  of  0.012  GPa,  for 
the  S+  and  S-  directions,  respectively. 

d)  Rotational  defects  in  RDX 

We  identified  a  new  type  of  point  defect  in  molecular  crystals:  rotational  defects.  In  this  case,  the  molecules  are  positioned  at 
the  proper  sites  of  the  crystal,  but  are  rotated  relative  to  their  neighbors.  A  rotational  defect  completely  breaks  the  point 
symmetry  at  the  respective  site,  which  is  different  from  other  types  of  molecular  rearrangement  reported  in  other  molecular 
crystals.  Such  rotated  molecules  are  stable  only  if  the  steric  hindrance  of  the  surrounding  molecules  helps  preserve  them  in  the 
rotated  positions. 

To  identify  whether  this  is  possible,  we  considered  a  perfect  RDX  crystal  in  which  a  molecule  was  rotated  to  various  positions, 
in  different  directions  relative  to  the  neighboring  molecules.  The  system  was  then  relaxed  to  determine  whether  the  rotated 
molecule  is  stable.  This  is  shown  schematically  in  Slide  13.  The  figure  represents  a  map  of  possible  molecular  rotations.  The 
rotations  are  characterized  by  the  angles  made  by  the  normal  to  the  ring,  as  shown  in  the  small  figure  next  to  the  map  on  slide 
13.  A  state  of  the  molecule  is  shown  by  a  point  on  this  map.  Upon  relaxation,  the  molecule  converts  to  another  structure  and 
orientation  -  this  process  is  represented  with  an  arrow  that  connects  the  initial  and  final  states  of  the  relaxation  process  in  this 
phase  space.  Note  that  several  attractors  appear,  i.e.  many  rotated  states  convert  to  a  small  number  (5)  of  states. 

These  rotated  states  are  stable  at  low  temperatures.  Upon  annealing  at  room  temperature,  all  these  rotational  defects 
disappear  and  the  perfect  crystal  is  recovered  (state  denoted  by  E).  We  computed  the  energy  barriers  for  the  formation  of  these 
defects  and  for  their  conversion  back  to  the  perfect  crystal  state.  The  barriers  for  defect  formation  are  on  the  order  of  1 .5  eV, 
while  those  for  their  annealing  (return  to  the  perfect  crystal  state)  are  about  0.3  eV.  The  formation  energy  of  the  various  defects 
is  in  the  vicinity  of  1  eV.  This  indicates  that:  a)  the  formation  energy  is  smaller  than  the  formation  energy  of  vacancies  (i.e.  the 
density  of  rotational  defects  should  be  larger  than  that  of  vacancies),  b)  thermal  fluctuations  are  insufficient  to  produce  rotational 
defects  at  all  temperatures  at  which  the  solid  phase  of  RDX  exists,  c)  the  0.3  eV  barrier  determines  a  characteristic  survival  time 
which  is  temperature  dependent  (i.e.  if  the  defect  is  formed  but  the  conditions  leading  to  its  formation  are  eliminated  at  t=0,  the 
time  interval  after  which  the  molecule  rotates  back  to  the  perfect  crystal  orientation  is  determined  by  this  “return  barrier”). 

Slide  14  shows  the  4  rotated  conformations  which  are  stable  in  addition  to  the  perfect  crystal  state.  These  correspond  to  the  5 
attractors  shown  in  the  map  of  slide  2  (labeled  A  to  E),  of  which  state  E  is  the  perfect  crystal  (not  shown  on  slide  14). 

Thermal  fluctuations  are  not  expected  to  produce  a  sizable  number  of  such  defects.  However,  we  have  identified  them  in 
dislocation  cores.  Two  dislocation  core  images  showing  rotated  molecules  (colored)  are  shown  on  slide  15.  Therefore,  we 
envision  that  rotational  defects  are  formed  in  cores  and,  as  dislocations  move  away,  are  left  behind  in  the  perfect  crystal. 
Thermal  fluctuations  anneal  them  back  in  a  characteristic  time  defined  by  the  0.3  eV  return  barrier;  this  time  is  on  the  order  of 
microseconds.  This  is  long  enough  for  these  defects  to  interact  with  phonons  and  other  dislocations  moving  in  the  same  glide 
plane.  The  importance  of  rotational  defects  to  plasticity  and  phonon  scattering  follows  from  this  observation. 

Relevance  for  initiation.  We  have  used  first  principle  calculations  (plane  wave  Density  Functional  Theory  with  projector 
augmented-wave  (PAW)  pseudopotentials  and  the  Perdew-Burke-Ernzerhof  (PBE)  exchange-correlation  functional 
implemented  within  the  generalized  gradient  approximation  (GGA)  in  VASP)  to  verify  the  atomistic  results  presented  before  and 
to  determine  the  electron  density. 

It  has  been  established  that  thermal  decomposition  in  solid  a-RDX  begins  with  N-N  bond  scission.  Therefore  we  examined  how 
electronic  charge  is  distributed  in  the  N-N02  region  of  the  defect  molecules.  If  charge  deficiency  over  the  N-N  bond  is 
considered  to  be  an  indicator  of  sensitivity,  one  needs  to  only  focus  on  the  most  charge  deficient  N-N  bond  in  the  molecule.  In 
addition,  since  only  one  N02  group  is  lost  per  molecule  during  decomposition,  the  N-N  bond  with  the  least  covalent  character 
can  be  expected  to  constitute  the  initial  fission  center.  Slide  16  shows  the  charge  distribution  along  2-D  slices  of  the  most 
deficient  N-N02  groups  for  the  defect  and  perfect  molecules.  For  the  perfect  crystal  case  (labelled  E),  the  charge  density  along 
an  axial  N-N02  group  is  shown,  and  one  can  see  that  the  electron  density  at  the  mid-point  of  the  N-N  bond  is  non-zero. 
However,  all  four  defects  (A-D)  show  some  charge  deficiency  at  the  mid-points  of  the  N-N02  bonds,  suggesting  that  these 
defect  molecules  may  be  more  sensitive  to  decomposition.  More  quantitatively,  in  the  perfect  crystal  (E)  the  minimum  density 
along  the  N-N  bond  is  0.327  e  Bohr-3.  For  the  defect  molecules  A-D,  the  N-N  electron  density  minimums  are  0.309,  0.288, 
0.318,  and  0.301  e  Bohr-3  respectively. 

How  to  identify  experimentally  rotational  defects?  To  assist  the  experimental  identification  of  rotational  defects,  we  computed 
the  eigenfrequencies  of  the  crystal  containing  the  4  types  of  rotational  defects  described  on  slide  14  in  a  region  of  high 
frequencies  where  the  spectrum  of  the  perfect  crystal  has  sparse  bands.  The  figure  on  slide  17  shows  the  perfect  crystal 
frequencies  with  black  bands  and  arrows.  The  colored  bands  correspond  to  the  various  defects.  As  seen,  band  splitting  occurs 
in  presence  of  rotational  defects  and  this  can  be  used  to  identify  them  using  Raman  spectroscopy. 

e)  Coarse  grained  models 

Two  models  of  the  RDX  (and  HMX)  molecule  exist  at  present  in  the  open  literature:  the  fully  explicit  atomistic  model  in  which  all 


atoms  are  represented  and  the  fully  rigid  molecule  developed  by  B.  Rice’s  group  at  ARL.  The  reference  model  is,  of  course,  the 
explicit,  atomistic  model.  It  has  been  shown  that  the  rigid  model  is  unable  to  capture  some  important  features  of  the  mechanics 
of  this  crystal.  Specifically,  the  predicted  stacking  fault  energies  are  incorrect  and  the  dislocation  cores  are  no  stabilized  by  this 
model.  On  the  other  hand,  any  coarse  grained  model  which  uses  larger  coarse  graining  levels  is  not  going  to  represent 
dislocations  explicitly  for  the  simple  reason  that  such  models  do  not  represent  explicitly  all  molecules  of  the  crystal.  Hence,  the 
coarse  graining  level  that  can  capture  the  atomistic  mechanics  right  and  hence  is  acceptable  should  be  lower  than  the  fully  rigid 
molecule. 

In  this  work  we  developed  a  family  of  models  that  interpolate  between  the  fully  flexible  (explicit  atomistic  model)  and  the  rigid 
molecule.  Their  ability  to  capture  specific  material  parameters  important  in  the  mechanics  of  the  crystal  has  been  tested.  Slide 
18  shows  the  5  models  of  this  family.  Model  I  corresponds  to  the  fully  flexible  molecule,  while  model  V  corresponds  to  the  fully 
rigid  molecule.  In  model  II  the  nitro  groups  are  rigid  and  the  bonds  connecting  the  nitro  group  to  the  ring  are  also  inextensible. 
The  rest  of  the  molecule  is  flexible.  Model  III  is  similar  to  model  II,  except  that  the  bonds  (C-C  and  C-N)  of  the  ring  are  also 
made  inextensible.  In  Model  IV,  the  whole  ring  is  made  rigid  while  the  rigid  nitro  groups  are  free  to  flop  about  the  ring.  In  Model 
V,  the  entire  molecule  is  a  rigid  object. 

Slide  19  shows  the  elastic  constants  of  the  crystal  predicted  with  the  5  models.  Models  II  and  III  give  results  close  to  the 
reference  (model  I).  Once  the  ring  is  made  rigid  (model  IV),  all  elastic  constant  depart  from  the  reference  by  about  10%.  The 
error  is  larger  when  the  rigid  model  is  used. 

The  lattice  parameters  are  less  sensitive  to  coarse  graining.  Models  III,  IV  and  V  lead  to  errors  (relative  to  the  reference,  model 

I)  of  about  2%.  The  pressure  dependence  of  the  elastic  constants  (relevant  for  shock)  is  the  same  for  all  models.  This  result  is 
shown  in  slide  20  for  one  of  the  lattice  parameters. 

The  coefficient  of  thermal  expansion  (CTE)  was  quite  sensitive  to  coarse  graining,  as  it  depends  on  phonon  mode  populations 
and  coarse-graining  directly  affects  the  normal  mode  distribution.  Models  II  and  III  gave  similar  values  (relative  to  the  reference) 
for  the  CTE  in  [100]  and  [010]  directions,  but  increased  by  about  20%  for  the  [001]  direction.  The  error  for  Models  IV  and  V  were 
more  pronounced,  going  as  low  as  50%  of  the  reference  CTE  value  for  the  [001]  direction.  Slide  21  shows  the  normalized  CTE 
for  the  5  models. 

These  models  indicate  qualitatively  the  degrees  of  freedom  that  need  to  be  active  for  accurate  thermal  and  elastic  behavior. 
Each  successive  model  offers  the  promise  of  lower  computational  cost 

(in  terms  of  degrees  of  freedom),  while  also  reducing  the  physics  that  can  be  captured  accurately.  In  this  context,  Model  III 
seems  to  offer  a  good  balance  of  features.  The  nitro  groups  can  be  coarsened  out  into  rigid  objects,  but  the  ring  has  to  stay 
flexible  (as  shown  by  poor  performance  of  models  IV  and  V).  This  study  provides  a  lower  bound  on  the  amount  of  coarse- 
graining  that  can  be  done  without  sacrificing  the  accuracy  with  which  the  physics  is  represented. 

Summary  statement 

To  summarize  the  main  advances  made  in  this  project  and  their  importance  in  multiscale  modeling: 

1 )  The  identification  of  active  slip  planes  using  atomistic  models  is  of  importance  for  crystal  plasticity  models  of  plastic 
deformation.  In  current  models,  plastic  deformation  is  considered  to  take  place  in  all  systems  of  the  respective  crystal  type.  We 
have  seen  that  on  some  systems  dislocations  are  simply  not  stable  and  hence  these  systems  should  not  be  considered  at  all. 

2)  Ranking  of  the  slip  systems  is  also  important  for  modeling  plastic  deformation  on  larger  scales.  Current  mesoscale  models 
have  no  input  regarding  the  ‘friction  stress’  (one  of  the  components  of  the  equation  of  evolution  of  plastic  strain  in  given  system) 
in  the  various  slip  systems.  This  work  provides  quantitative  information  in  this  sense. 

3)  Rotational  defects  were  described  for  the  first  time.  We  quantified  their  potential  importance  in  initiation. 

4)  A  family  of  coarse  grained  models  has  been  developed  for  RDX  and  it  was  shown  that  coarse  graining  the  entire  molecules 
in  blobs  does  not  lead  to  accurate  results.  As  long  as  discrete  coarse  grained  models  are  to  be  used,  model  III  described  above 
is  the  highest  level  of  coarse  graining  acceptable. 


Technology  Transfer 

This  work  was  performed  in  collaboration  with  the  groups  of  Dr.  Betsy  Rice  and  Peter  W.  Chung  from  the  Army  Research  Labs. 
Dr.  Rice  is  specialized  in  quantum  mechanical  and  atomistic  modeling  of  energetic  materials.  Our  contribution  was  to  add 
physical  understanding  and  to  extend  the  range  of  scales  accessible  to  their  models,  which  include  atomistic  (nano)  and  coarse 
grained/mesoscale  (micro)  models.  We  interacted  closely  with  Dr.  Peter  Chung  (from  the  Computational  &  Information  Sciences 
Directorate  of  ARL  until  Spring  2014)  who  has  contributed  significantly  to  the  work  on  the  Peierls  stress  of  dislocations  in  RDX 
which  has  been  published  in  201 3.  This  work  initiated  as  a  result  of  our  visit  to  ARL  in  April  2011  and  was  discussed  in 
subsequent  meetings  at  RPI.  Dr.  Chung  has  served  in  the  PhD  committee  of  Nithin  Mathew.  Nithin  graduated  in  January  2013. 

In  April  2014  we  visited  the  group  of  Dr.  Betsy  Rice  at  ARL  and  gave  a  talk  on  the  results  of  the  work  reported  here.  This  led  to 
the  invitation  to  submit  a  proposal  for  the  initiation  of  a  Cooperative  Agreement  with  ARL.  The  proposal  has  been  approved.  The 
Cooperative  Agreement  provides  a  framework  for  sharing  information  (model  development  and  physical  results)  between  ARL 
and  RPI  and  provides  access  to  the  RPI  group  to  the  ARL  supercomputing. 

We  have  also  interacted  with  Dr.  Plimpton  from  Sandia  National  Labs.  Dr.  Plimption  is  the  primary  developer  of  the  LAMMPS 
molecular  dynamics  software.  We  worked  together  on  specific  extensions  of  LAMMPS. 

Both  the  PI  and  the  Co-Pi  contributed  in  the  past  to  another  project  (STTR)  involving  ARL  and  Simmetrix  Inc.,  project  focused 
on  the  development  of  multiscale  modeling  technologies  for  nanomaterials.  The  tools  developed  part  of  this  project  are  currently 
implemented  at  ARL  and  have  been  used  in  the  current  work  both  at  ARL  and  RPI. 

The  two  Pis  worked  jointly  on  a  project  sponsored  by  IBM  focused  on  the  development  and  use  of  multiscale  technologies  for 
the  purpose  of  identifying  conditions  for  defect  nucleation  in  microelectronic  devices.  The  broader  goal  of  this  set  of 
technologies  is  to  bridge  the  gap  between  fabrication  and  device  performance  and  reliability  (“front  end”  and  “back  end”). 

The  Pis  were  also  involved  in  the  past  with  two  NSF-funded  projects  on  multiscale  modeling  of  polymer-based  nanocomposites, 
each  funded  over  a  duration  of  3  years  at  levels  above  $500,000. 

A  three-year  NSF  project  ($325,000  total  funding)  on  the  continued  development  of  the  RPI’s  AMSI  multiscale  computing 
framework  ends  in  July  2015. 


Energetic  materials  are  molecular  crystals 


Chemical  structure  of  the  RDX  molecule  Structure  of  a-RDX  crystal 


Single  molecule 
Red  =  Oxygen 
Blue  =  Nitrogen 
Cyan  =  Carbon 
White  =  Hydrogen 


c) 


(orthorombic  -  [100]  projection) 


(a)  Conformational  defects  in  RDX 
-what  are  conformers?- 


•  Estimation  of  energy  barriers 


Caee  conformer 


(a)  Conformational  defects  in  RDX 
-in  a-RDX  a  fraction  of  molecules  are  conformers  - 


•  a-RDX  is  composed  from  molecules  in  the  Caae 
state. 

•  Under  uniaxial  tensile  (and  compressive)  loading 
in  the  [001]  direction  (as  indicated  in  the  figure),  a 
fraction  of  Caae  molecules  convert  to  Caee. 


Concentration  of  Caee  after  imposing  a  step  strain  of 
0. 7%  in  [001].  The  inset  shows  the  dependence  of  the 
stable  Caee  concentration  (fraction)  with  strain. 
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(a)  Conformational  defects  in  RDX 
-conformers  tend  to  cluster  in  space  and  transitions 

tend  to  cluster  in  time  - 


Transformation  from  Caae  to  Caee  and  back  is  a  stochastic  process.  However, 
the  events  are  spatially  and  temporally  correlated.  This  is  expected  to  be 
important  for  the  mechanism  of  thermal  energy  concentration  and  hot  spot 
formation. 


Pair  distribution  function  of  Caee  molecules 
indicating  clustering  of  transformation  events  in 
the  crystal. 


Correlation  integral  analysis  of  time  series  of 
Caae  to  Caee  transitions  indicating  temporal 
correlations.  The  range  of  these  temporal 
correlations  is  approximately  1  ps. 


(a)  Conformational  defects  in  RDX 
-a  transition  leads  to  local  temperature  fluctuations- 


“Temperature”  of  intramolecular  vibration  modes  following 
conformation  transitions  of  a  single  molecule  in  the  RDX 

crystal 


t(ps) 


(b)  Stability  of  dislocations  in  RDX 
-family  of  crystallographic  planes  on  which  slip  occurs- 


+  001 
_  0  10 


•Slip  on  (001),  {021}, 
{011}  planes  has  been 
reported  in  experimental 
studies.  Note  that  direct 
observations  of 
dislocations  in  RDX  is  not 
possible.  Computer 
simulations  can  clarify  the 
preferred  glide  systems, 
the  mobility  in  various 
glide  systems,  the  Peierls 
barriers  to  motion. 


(100)  projection  of  RDX  crystal  structure 
indicating  the  slip  planes  in  which  glide  can 
happen  in  the  [100]  direction  (perpendicular  to 
the  figure). 
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(b)  Stability  of  dislocations  in  RDX 
-specific  stable  types  of  dislocations- 


Structure  of  the  (010) [100]  edge  dislocation,  (a)  The  dislocation  with  b  =  [100] 
splits  into  two  partials  of  b  =  1/2[100].  (b)  Structure  of  the  b  =  1/2[100]  partial. 


The  (010)[100]  screw  configuration  was  found  to  be  unstable. 
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(b)  Stability  of  dislocations  in  RDX 
-specific  stable  types  of  dislocations- 


i  1  i  1  - T 


Core  structure  of  the  (001) [010]  edge  dislocation  with  b  =  1/2[010]. 

(a)  Relative  displacement  of  COM  across  the  glide  plane.  Note  the  occurrence  of  the 
screw  component  (blue  curve) 

(b)  Molecular  conformational  state  of  the  core  in  (a).  Caaa  molecules  are  colored  in  red, 

Caee  are  in  blue,  and  Caae  are  in  gray.  8 


u(x)/l100 


(b)  Stability  of  dislocations  in  RDX 
-specific  stable  types  of  dislocations- 


— o—  Volterra  field 

I 

• 

—a—  Volterra  field 

—•—Atomistic  solution  [100] 

0.6- 

—•—Atomistic  solution  [100] 

—o— Atomistic  solution  [0-12] 

— o—  Atomistic  solution  [0-12] 

x/l012  x/l100 


Core  structure  of  b  =  1/2  [100]  dislocations  residing  on  the  (021)  plane, 
(a)  Screw  dislocation  (b)  Edge  dislocation 


Peierls  Stress  (GPa) 


(c)  Ranking  of  slip  systems  in  RDX 
-Peierls  stress  values- 


Labels  “e”  and  “s”  stand  for  edge 
and  screw  dislocations.  The  bars 
represent  the  range  in  which  the 
Peierls  stress  is  observed.  The 
limits  defining  the  two  shaded 
areas  are  selected  arbitrarily  to 
correspond  to  Gv/20  and  Gv/10, 
where  Gv  is  the  Voigt  shear 
modulus. 
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(c)  Ranking  of  slip  systems  in  RDX 
-Implications  of  the  ranking- 


Activation  of  slip  systems  for  (111) 
and  (210)  shocks 


Activation  of  slip  systems  in  indentation 
experiments  perpendicular  to  the  (001),  (021), 
and  (210)  crystal  planes 


(c)  Ranking  of  slip  systems  in  RDX 
-Slip  asymmetry- 


(d)  Rotational  defects  in  RDX 
Identification  of  stable  rotated  molecules 


180  - 1 - 1 - 1 - I - r 
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Map  of  the  conformation  (rotation)  space.  The  arrows  show  how  a  molecule  in  a 
certain  initial  rotation  state  relaxes.  The  points  to  which  multiple  arrows  point  are 
attractors.  All  states  labeled  A  to  E  are  stable  at  low  temperature,  but  convert  to  the 
perfect  crystal  orientation  upon  annealing  at  room  temperature. 


(d)  Rotational  defects  in  RDX 
-Four  stable  rotational  defects- 


Four  configurations  of  rotated  molecules  (rotational  point  defects)  which  are  stable  at 
low  temperatures  and  which  correspond  to  the  4  attractors  of  the  map  of  slide  13.  The 
fifth  point  on  that  map  represents  the  perfect  crystal  orientation. 
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(d)  Rotational  defects  in  RDX 
-Rotational  defects  are  present  at  dislocation  cores- 


Edge  and  screw  dislocation  cores  and  rotational  defects  shown  in  color  in  the  core 
region.  These  defects  are  stable  in  the  cores  even  at  room  temperature.  When  the 
dislocation  moves  away,  the  defects  heal  back  and  disappear.  This  process  takes  a 
time  comparable  with  the  time  required  for  the  next  dislocation  moving  in  the  same 
glide  plane  to  arrive  at  the  same  location.  Hence,  we  conjecture  that  these  structures 
interact  strongly  with  dislocations.  i 5 


?n  (d)  Rotational  defects  in  RDX 

'ectron  charge  density  redistribution  in  rotational 

defects- 


Electron  charge  distribution  across  the  most  charge  deficient  A/-A/O2  groups 
in  the  perfect  crystal  (E)  and  Defect  (A-D)  molecules.  The  color  scales  from  0  to  1.32  (Bohr3)  with 
contours  at  0.45  to  1.20  in  multiples  ofO.  15.  All  defects  display  density  reduction  at  the  mid-point  of  the 

N-N  bond,  with  Defect  B  having  the  largest  deficiency. 


(d)  Rotational  defects  in  RDX 
■Identification  of  rotational  defects  using  Raman ■ 
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Splitting  of  high  frequency  Raman  active  modes  in  the  presence  of  rotational  defects  of 
type  A  to  D.  The  modes  of  defective  crystals  are  shown  in  color,  while  those  of  the 
perfect  crystal  are  shown  in  black  (and  arrows) 


(e)  Coarse  grained  models  of  RDX  crystal 
-Hierarchy  of  models  considered- 


Model  I:  Fully  flexible 


Model  II:  Ring  flexible,  N-N 
bonds  rigid  and  nitro  groups  rigid 


Increasing  coarse  graining  level 


Model  III:  Ring  flexible,  but  bonds  Model  IV:  Everything  rigid  except  Model  V-  Fully  rigid 
very  stiff  axially.  Nitro  groups  rigid  the  nitro  groups  (which  are  rigid) 
but  can  move  relative  to  the  ring.  can  move  relative  to  the  ring 


Increasing  coarse  graining  level 


Hierarchy  of  coarse  grained  models  developed  in  this  work 


(e)  Coarse  grained  models  of  RDX  crystal 
-Effect  of  coarse  graining  on  elastic  constants- 
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Prediction  of  crystal  elasticity  using  the  6  models  defined  in 

slide  5. 


19 


Lattice  Constant  a  (angstrom) 


(e)  Coarse  grained  models  of  RDX  crystal 
-Effect  of  coarse  graining  on  lattice  parameters  and 

bulk  modulus- 


Lattice  constant  versus  pressure 
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Variation  of  one  of  the  lattice  parameters  (in  Angstroms)  with 
pressure,  as  predicted  by  the  5  coarse  grained  models. 
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(e)  Coarse  grained  models  of  RDX  crystal 
-Effect  of  coarse  graining  on  thermal  expansion- 


Normalized  CTE  for  250-300K 
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■  I  ■  II  ■  III  ■  IV  ■  V  ■ 

Variation  of  the  CTE  (normalized  with  reference  model  I)  in 
three  principal  directions  of  the  crystal  at  zero  pressure  and  for 
all  coarse  grained  models. 
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